This file contains models examining the functional relationship between agricultural landscape diversity and agricultural yields. This work applies hierachical Bayesian spatiotemporal modeling techniques to estimation of the effect of diversity of the agricultural landscape on the yield of corn, soy, and winter wheat.

REPRESENTING SPACE

The yield data available is at a county scale and the distribution of yields across space exhibits strong autocorrelation where yields in neighboring counties are more alike than yields in distant counties. This spatial autocorrelation is accounted for using a standard Conditional Autoreggressive dependency model based on adjacency for all counties in the conterminous US. In order to account for additional county-specific factors that contribute to yields a county iid random effect term is also included, yielding a Besag-York-Mollie (BYM) spatial dependency model. A county adjacency matrix is created for each crop under investigation. For regional models a seperate county adjacency matrix for each region-crop combination is created.

MODELS

The model form is a random-effects panel model with county spatial effects, non-parametric SDI, non-parametric Temp and Precip controls, and a quadratic time trend (similar to Schlenker and Roberts, 2009) for each agro-eco-region. The spatial effects account for standard, unstructured or iid, county random effects (deviation of the county mean from the population mean) as well as conditional autogregressive (CAR) structured residuals between counties. These effects account for time-invariant factors that differ across counties (however assume that the effects are not correlated with errors hence do not fully account for potential ommitted variables). The quadratic time trends for each agro-eco-region account for factors that are space-invariant (within an agro-eco-region) and vary across time (without the overfitting that was observed with the use of the rw2 random time effects). This effect is intended to account for factors that produce changes over time, but that may differ across agro-eco-regions. These effects are expected to partially account for differences in technology use and adoption rates over time in different regions and differences in agricultural management styles and techniques across regions that result in temporally varying responses to climate conditions. Note that while there clearly are differences in management styles, climate responses, and technology use within agro-eco-regions these effects account for aspects of the agro-eco-region as a whole that constrain or enable certain types of actions.

formula<-YIELD ~ 1 + f(TP, model=“rw1”,scale.model=T) + f(SDD,model=“rw1”,scale.model=T) + f(GDD,model=“rw1”, scale.model=T) + f(SDI, model=“rw1”, scale.model=T)+ PERC_IRR + ACRES + f(CNTY, model=“bym”, graph=H.adj.corn, scale.model=TRUE) + f(AERCODE.id, YEAR.id)+ f(AERCODE.id2, YEAR.id2)

Previous models employed r-inla default prior settings for the random effects and used a reduced precision fixed effect prior. This round of models keeps the reduced precision fixed effect prior and further employs penalized complexity priors (https://arxiv.org/pdf/1403.4630v4.pdf) that use a scaling factor to specify priors based on the limits of the data creating a weakly informative prior. Default and recomended settings for PC priors as provided in the R-INLA documentation and in https://arxiv.org/pdf/1403.4630v4.pdf are employed. These models also employ sum-to-zero constraints for all random effects(https://link.springer.com/article/10.1007%2Fs00477-017-1405-0).

ABOUT DIAGNOSTICS

In the returned model ouputs several diagnostics are reported. Some of these diagnostics are computed by R-INLA, and additional posterior predicve and cross-validation checks were manually computed. Descriptions of the diagnostics are adapted from Spatial and Spatio-temporal Bayesian Models with R-INLA (2015) by Blangiardo and Cameletti.

Manually calculated posterior predictive checks (PPC) were obtained using the predictive distribution estimated by R-INLA. These checks use all observations for both model evaluation and predictive checks (hence are not a form of cross-validation).

In addition, an independent cross-validation model was run using 6/7ths of the data for model fit estimation and reserving a held-out sample of 1/7th of the data for comparison. Three diagnostic metrics (MSE, R2, and NSE) were calculated for the cross-validation models. As a smaller sample size is used for these caluculations they are more sensitive to outliers than the previously reported diagnostics. Therefore, outlier error values greater than twice the mean of the observed data were removed prior to calculation of these diagnostics and the count of outliers is also reported.

US EXTENTS CORN MODELS

This section estimates the model described above for corn for three metrics of agricultural landscape diversity (SIDI, SDI, and RICH). Model estimates, diagnostics, functional response plots, and maps of spatial effects are all reported.

[[1]]


Table: Summary Table of Model Estimates

                                                mean          sd   0.025quant    0.5quant   0.975quant
----------------------------------------  ----------  ----------  -----------  ----------  -----------
(Intercept)                                   4.0998      0.0218       4.0569      4.0998       4.1425
PERC_IRR                                      0.0091      0.0006       0.0079      0.0091       0.0102
ACRES                                         0.0000      0.0000       0.0000      0.0000       0.0000
Precision for the Gaussian observations      29.9498      0.4522      29.0791     29.9435      30.8537
Precision for TP                            323.2940    153.4361     110.6326    295.8486     700.5702
Precision for SDD                            20.8065      5.5533      11.9756     20.1122      33.6745
Precision for GDD                           114.8603     46.6638      50.5940    105.9077     230.4463
Precision for SIDI                         2003.0531   1482.3413     440.7825   1609.4662    5900.7110
Precision for CNTY                           22.2986      1.1455      20.1984     22.2404      24.7028
Phi for CNTY                                  0.9990      0.0017       0.9941      0.9997       1.0000
Precision for AERCODE.id                    355.2637     81.5375     221.4464    346.3688     540.6241
Precision for AERCODE.id2                  9911.1632   2139.6328    6373.4412   9685.9031   14747.1260

[[2]]


[[3]]


Table: Model Diagnostic Metrics

      DIC        CPO         MSE          R2
---------  ---------  ----------  ----------
 -4756.95   2209.102   0.0289812   0.7082578
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.7857263 0.0382403 0.6607148 2
[[1]]


Table: Summary Table of Model Estimates

                                                mean          sd   0.025quant    0.5quant   0.975quant
----------------------------------------  ----------  ----------  -----------  ----------  -----------
(Intercept)                                   4.1323      0.0212       4.0908      4.1323       4.1740
PERC_IRR                                      0.0088      0.0006       0.0076      0.0088       0.0100
ACRES                                         0.0000      0.0000       0.0000      0.0000       0.0000
Precision for the Gaussian observations      29.9944      0.4531      29.1064     29.9934      30.8910
Precision for TP                            323.9264    157.6665     113.1573    292.9302     715.6133
Precision for SDD                            20.6454      5.5220      11.8783     19.9501      33.4501
Precision for GDD                           112.2088     45.1103      49.0859    103.9458     223.3186
Precision for SDI                           885.7269    600.4967     230.4002    729.7607    2465.8757
Precision for CNTY                           22.0985      1.1524      19.7891     22.1258      24.3032
Phi for CNTY                                  0.9993      0.0013       0.9955      0.9999       1.0000
Precision for AERCODE.id                    371.9377     85.8291     231.0064    362.6058     567.0494
Precision for AERCODE.id2                  9913.3866   2181.6194    6380.7361   9652.6657   14916.9165

[[2]]


[[3]]


Table: Model Diagnostic Metrics

       DIC        CPO         MSE          R2
----------  ---------  ----------  ----------
 -4772.191   2214.884   0.0289038   0.7084921
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.7397979 0.0423634 0.6614034 2
[[1]]


Table: Summary Table of Model Estimates

                                                 mean          sd   0.025quant     0.5quant   0.975quant
----------------------------------------  -----------  ----------  -----------  -----------  -----------
(Intercept)                                    4.1761      0.0230       4.1310       4.1760       4.2213
PERC_IRR                                       0.0084      0.0006       0.0073       0.0084       0.0096
ACRES                                          0.0000      0.0000       0.0000       0.0000       0.0000
Precision for the Gaussian observations       29.9331      0.4532      29.0635      29.9248      30.8437
Precision for TP                             341.6731    174.5707     118.6400     304.0864     785.2311
Precision for SDD                             21.8095      5.7672      12.5552      21.1180      35.1047
Precision for GDD                            117.2635     47.1568      51.7564     108.4094     233.3232
Precision for RICH                          1079.5893    566.2113     350.2677     960.4435    2510.8869
Precision for CNTY                            22.4371      1.1727      20.0898      22.4638      24.6855
Phi for CNTY                                   0.9992      0.0014       0.9951       0.9999       1.0000
Precision for AERCODE.id                     412.9991     97.1548     248.1269     404.8015     628.0932
Precision for AERCODE.id2                  11123.1604   2410.7423    6922.2920   10965.9996   16340.8399

[[2]]


[[3]]


Table: Model Diagnostic Metrics

       DIC        CPO         MSE          R2
----------  ---------  ----------  ----------
 -4762.183   2213.329   0.0289657   0.7080809
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.7112438 0.0413223 0.6744548 2

US EXTENTS SOY MODELS

This section estimates the model described above for soy for three metrics of agricultural landscape diversity (SIDI, SDI, and RICH). Model estimates, diagnostics, functional response plots, and maps of spatial effects are all reported.

[[1]]


Table: Summary Table of Model Estimates

                                                 mean           sd   0.025quant    0.5quant   0.975quant
----------------------------------------  -----------  -----------  -----------  ----------  -----------
(Intercept)                                    3.2781       0.0141       3.2504      3.2781       3.3057
PERC_IRR                                       0.0065       0.0005       0.0054      0.0065       0.0075
ACRES                                          0.0000       0.0000       0.0000      0.0000       0.0000
Precision for the Gaussian observations       49.7800       0.7936      48.2088     49.7863      51.3274
Precision for TP                             695.5598     493.8363     133.6313    574.7619    1977.3435
Precision for SDD                             72.6201      21.8427      38.2643     69.8419     123.2019
Precision for GDD                            287.6908     106.0468     123.2948    274.2293     533.0373
Precision for SIDI                         15200.4452   22564.7069    1532.3233   8593.0572   69739.3076
Precision for CNTY                            33.1739       1.7864      29.8960     33.0843      36.9240
Phi for CNTY                                   0.9974       0.0031       0.9889      0.9986       1.0000
Precision for AERCODE.id                     250.6527      55.3843     158.1464    245.1963     375.1507
Precision for AERCODE.id2                   7275.8801    1535.8317    4640.2622   7154.3799   10643.1035

[[2]]


[[3]]


Table: Model Diagnostic Metrics

       DIC        CPO        MSE         R2
----------  ---------  ---------  ---------
 -9225.362   4393.392   0.017492   0.751808
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.7736001 0.0258817 0.7081603 2
[[1]]


Table: Summary Table of Model Estimates

                                                mean          sd   0.025quant    0.5quant   0.975quant
----------------------------------------  ----------  ----------  -----------  ----------  -----------
(Intercept)                                   3.2882      0.0131       3.2624      3.2882       3.3138
PERC_IRR                                      0.0064      0.0005       0.0054      0.0064       0.0074
ACRES                                         0.0000      0.0000       0.0000      0.0000       0.0000
Precision for the Gaussian observations      49.8351      0.8359      48.1051     49.8774      51.3696
Precision for TP                            663.5682    465.6018     123.1728    552.2777    1868.7598
Precision for SDD                            72.3292     21.7891      37.3316     69.8613     122.2682
Precision for GDD                           285.3265    105.9807     121.5044    271.7122     530.9557
Precision for SDI                          2760.1632   1710.5156     626.3132   2390.6095    7070.1910
Precision for CNTY                           33.5448      2.2422      28.8027     33.7247      37.4636
Phi for CNTY                                  0.9973      0.0036       0.9872      0.9988       1.0000
Precision for AERCODE.id                    259.5931     61.5184     165.0557    250.4716     404.9463
Precision for AERCODE.id2                  7220.0371   1633.9353    4664.8174   6992.7377   11034.4481

[[2]]


[[3]]


Table: Model Diagnostic Metrics

      DIC        CPO         MSE        R2
---------  ---------  ----------  --------
 -9241.69   4401.045   0.0174666   0.75212
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.733186 0.0242373 0.7366562 2
[[1]]


Table: Summary Table of Model Estimates

                                                mean          sd   0.025quant    0.5quant   0.975quant
----------------------------------------  ----------  ----------  -----------  ----------  -----------
(Intercept)                                   3.2933      0.0133       3.2672      3.2933       3.3194
PERC_IRR                                      0.0064      0.0005       0.0054      0.0064       0.0074
ACRES                                         0.0000      0.0000       0.0000      0.0000       0.0000
Precision for the Gaussian observations      49.8778      0.7920      48.3326     49.8735      51.4507
Precision for TP                            679.8154    524.4508     139.9951    538.2348    2061.3030
Precision for SDD                            71.6777     21.6642      37.6102     68.9289     121.7865
Precision for GDD                           284.4466    105.1302     125.8213    269.2925     533.5239
Precision for RICH                         4902.6156   3740.9666    1252.4640   3843.1526   14820.7651
Precision for CNTY                           33.0099      1.7383      29.7691     32.9423      36.6026
Phi for CNTY                                  0.9981      0.0024       0.9914      0.9991       1.0000
Precision for AERCODE.id                    245.0916     56.0029     155.4249    238.0331     374.5956
Precision for AERCODE.id2                  7383.8038   1557.0342    4647.7376   7289.4864   10732.7118

[[2]]


[[3]]


Table: Model Diagnostic Metrics

       DIC        CPO         MSE          R2
----------  ---------  ----------  ----------
 -9242.452   4399.599   0.0174389   0.7522154
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.7599772 0.0246851 0.7117944 2

US EXTENTS WINTER WHEAT MODELS

This section estimates the model described above for winter wheat for three metrics of agricultural landscape diversity (SIDI, SDI, and RICH). Model estimates, diagnostics, functional response plots, and maps of spatial effects are all reported.

[[1]]


Table: Summary Table of Model Estimates

                                                mean         sd   0.025quant    0.5quant   0.975quant
----------------------------------------  ----------  ---------  -----------  ----------  -----------
(Intercept)                                   3.7795     0.0297       3.7217      3.7793       3.8384
PERC_IRR                                      0.0069     0.0006       0.0057      0.0069       0.0081
ACRES                                         0.0000     0.0000       0.0000      0.0000       0.0000
Precision for the Gaussian observations      30.3264     0.5634      29.1978     30.3369      31.4122
Precision for TP                            194.1723    82.3368      74.5400    181.1520     392.3227
Precision for SDD                           659.1264   310.7045     257.3424    592.3727    1448.5368
Precision for GDD                           531.0536   277.8829     178.1044    470.7707    1239.1956
Precision for SIDI                          926.9797   648.8700     229.7683    756.4284    2633.5675
Precision for CNTY                           12.4511     0.8780      10.7602     12.4427      14.2149
Phi for CNTY                                  0.9756     0.0109       0.9499      0.9772       0.9919
Precision for AERCODE.id                    134.2762    26.6501      89.7140    131.6345     194.0044
Precision for AERCODE.id2                  4486.9209   898.7138    3039.6699   4374.9881    6548.7653

[[2]]


[[3]]


Table: Model Diagnostic Metrics

       DIC        CPO         MSE          R2
----------  ---------  ----------  ----------
 -3410.903   1585.825   0.0281662   0.7687149
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.7669969 0.0404798 0.7480352 2
[[1]]


Table: Summary Table of Model Estimates

                                                mean         sd   0.025quant    0.5quant   0.975quant
----------------------------------------  ----------  ---------  -----------  ----------  -----------
(Intercept)                                   3.8526     0.0317       3.7911      3.8523       3.9153
PERC_IRR                                      0.0063     0.0006       0.0051      0.0063       0.0075
ACRES                                         0.0000     0.0000       0.0000      0.0000       0.0000
Precision for the Gaussian observations      30.2602     0.5574      29.1759     30.2563      31.3698
Precision for TP                            186.7372    82.9091      75.7365    170.0672     393.9371
Precision for SDD                           615.3965   278.2376     241.1038    560.4547    1310.6071
Precision for GDD                           501.0939   254.4968     174.8109    446.5999    1147.2450
Precision for SDI                           403.7058   324.3212      90.2708    312.1551    1261.6378
Precision for CNTY                           12.5667     0.8843      10.9004     12.5447      14.3715
Phi for CNTY                                  0.9750     0.0110       0.9488      0.9766       0.9912
Precision for AERCODE.id                    140.7010    27.9104      93.6204    138.0694     203.0532
Precision for AERCODE.id2                  4584.8992   906.0762    3101.9769   4480.3011    6653.0673

[[2]]


[[3]]


Table: Model Diagnostic Metrics

       DIC        CPO         MSE          R2
----------  ---------  ----------  ----------
 -3397.041   1578.491   0.0281994   0.7682754
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.7362974 0.0391073 0.7371566 2
[[1]]


Table: Summary Table of Model Estimates

                                                mean         sd   0.025quant    0.5quant   0.975quant
----------------------------------------  ----------  ---------  -----------  ----------  -----------
(Intercept)                                   3.8501     0.0300       3.7917      3.8499       3.9095
PERC_IRR                                      0.0063     0.0006       0.0051      0.0063       0.0075
ACRES                                         0.0000     0.0000       0.0000      0.0000       0.0000
Precision for the Gaussian observations      30.4295     0.5616      29.3509     30.4205      31.5553
Precision for TP                            198.8608    88.1498      79.5291    181.6452     418.3981
Precision for SDD                           660.4567   308.0338     256.7287    596.0708    1436.6173
Precision for GDD                           530.9166   290.7141     181.0598    461.9682    1283.2551
Precision for RICH                          568.7480   257.2415     218.6604    519.3383    1208.6633
Precision for CNTY                           12.6909     0.9147      10.9541     12.6714      14.5458
Phi for CNTY                                  0.9731     0.0116       0.9459      0.9748       0.9904
Precision for AERCODE.id                    136.0891    26.7152      90.2826    133.8813     194.9322
Precision for AERCODE.id2                  4616.9970   879.1737    3107.1747   4544.2102    6555.3250

[[2]]


[[3]]


Table: Model Diagnostic Metrics

       DIC        CPO         MSE          R2
----------  ---------  ----------  ----------
 -3437.564   1598.147   0.0280755   0.7692436
[[1]]


[[2]]

[[1]]

[[1]]

[[1]]

[[1]]

Cross Validation Diagnostic Metrics
R2 MSE NSE OutlierCount
0.7940621 0.0404092 0.7583557 2

SUMMARY PLOTS

Summary plots showing the effect of the smooth, random-walk, effects of diversity and climate on log(Yield) for each model. All crops are presented on a single plot, and new plots are created for different diversity metrics and different climate variables. Solid lines represent the median effect and the colored bands represent the 95% credibility limits as obtained from the estimated posterior distribution.

Effect of Diversity Metrics

Effect of Climate Variables